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CONVERGENCE ESTIMATES FOR MULTIDISCIPLINARY ANALYSIS AND 

OPTIMIZATION 

EYAL A RIAN * 


Abstract. A quantitative analysis of coupling between systems of equations is introduced. This analysis 
is then applied to problems in multidisciplinary analysis, sensitivity, and optimization. For the sensitivity 
and optimization problems both multidisciplinary and single discipline feasibility schemes are considered. 
In all these cases a “convergence factor” is estimated in terms of the Jacobians and Hessians of the system, 
thus it can also be approximated by existing disciplinary analysis and optimization codes. The convergence 
factor is identified with the measure for the “coupling” between the disciplines in the system. Applications 
to algorithm development are discussed. Demonstration of the convergence estimates and numerical results 
are given for a system composed of two non-linear algebraic equations, and for a system composed of two 
PDEs modeling aeroelasticity. 

Key words, convergence, MDO, sequential, Gauss-Seidel, adjoint, Hessian 

Subject classification. Applied Numerical Mathematics 

1. Introduction. In the last three years there have been a growing interest in the engineering com- 
munity in the numerical solution of “multidisciplinary” (MD) problems [l]-[3], A MD problem is defined 
in this work to be a problem which involves the solution of a set of equations which can be divided into at 
least two sets (disciplines). Three different MD problems are considered: the “Multidisciplinary Analysis” 
problem (MDA), the “Multidisciplinary Sensitivity” (Adjoint) problem (MDS), and the “Multidisciplinary 
Optimization” problem (MDO). The MDA problem is to solve a system of state equations, the MDS problem 
is to compute the sensitivity of some given quantity with respect to “design” parameters which appear in 
the state equations, and the MDO problem is to minimize a given cost- functional with respect to the design 
variables. Typically each set of state equations models a different physical process and therefore many times 
each set is dominated by different scales; the Jacobian of each set of equations is characterized by a differ- 
ent eigenvalue spectra resulting in a large and ill-conditioned problem. Also the MD problem introduces a 
practical difficulty of connecting large codes of different origin. Whenever possible (computational- wise) it 
is desirable to perform the computation in a “loosely coupled” scheme in which most of the computation is 
done in the disciplinary levels. 

The recent effort in MDO is aiming at the development of new strategies and algorithms to solve the 
MD problems in a “tightly- coupled” approach. In some cases the optimization process is performed in the 
“system level” while the analysis and sensitivity analysis are performed in the disciplinary levels. Having 
the sensitivities of all the disciplines, the system optimization iteration convergence can be analyzed with 
classical methods developed for single discipline optimization. For such an approach this work is relevant 
only to the development of algorithms to solve the MDA and MDS problems. References [l]-[3] summarize 
most of the work done so far in MDO. R.J. Balling and J. Sobieski [4] propose several multi-level algorithms 
for the solution of MDO problems in which the design process is decomposed into the disciplinary levels. 
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VA 23681-0001 (email: arian@icase.edu). This research was supported by the National Aeronautics and Space Administration 
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E.J. Cramer et. al. [5] propose three fundamental formulations to MDO based on feasibility requirements: 
Multidisciplinary Feasible (MDF), Individual Discipline Feasible (IDF), and All-At-Once (A AO). These 
formulations were applied by G.R. Shubin [6] to aeroelasticity and it is reported that in general the less 
feasible formulation, A AO, was the most efficient in practice. 

The classical approach to such problems is the sequential ( “one- at- a- time” ) approach; solve each disci- 
plinary problem separately, while freezing the variables of the other disciplines. Up-to-date the sequential 
approach is used in applications to solve multidisciplinary optimization problems, for example in aircraft 
design (e.g.,[7j). The sequential approach has the advantage that it is the simplest to apply since typically 
there is a great deal of information and experience on numerical methods for the solution of each of the 
disciplinary problems. However, in many cases such an approach is undesirable since it might not converge 
to the MD solution at all [8], and in some other cases it might converge to the solution but slowly, requiring 
numerous disciplinary solutions and exchange of data between the disciplines at each MD iteration. The 
sequential approach is known in iterative methods for solving algebraic equations as the Gauss- Seidel algo- 
rithm (for example [9] p.508), i.e., each of the equations in the system is solved separately for one variable 
while freezing the variables of the other equations. Therefore, in this paper such an approach will be called 
“a Generalized Gauss-Seidel” (GGS) approach. The GGS approach has a parallel version in which the disci- 
plinary problems are solved separately in parallel rather than sequentially. This scheme is known in iterative 
methods for solving algebraic equations as the Jacobi algorithm and therefore in this paper such an approach 
will be called the “Generalized Jacobi” (GJ) approach. 

Qualitative predictions of the GGS convergence for an MDO problem were done in [10], where it was 
shown that the MD Hessian determines the coupling in a static aeroelastic optimization problem. Local 
mode analysis of the MD Hessian was performed to prove that for the non-smooth error mode components 
the problem has very weak coupling between aerodynamic design and structural design and therefore it 
predicts that a sequential approach will be efficient for that problem. 

In this paper, a quantitative convergence analysis is developed for the GGS (sequential), and GJ (par- 
allel), schemes applied to MDA, MDS, and MDO problems. An upper bound for the convergence rate of 
the GGS scheme is estimated by the maximal singular value of a “convergence matrix”, .M, that is given in 
terms of the Jacobians and Hessians of the problem. In problems governed by PDEs the convergence matrix 
is an operator (typically non-differential) which determines the “coupling” between the disciplines in the 
system. We claim that these operators have a fundamental role in the mathematical theory of MDO and 
we give an explicit analysis of these operators for an aeroelastic model problem (Sec. 6). The motivation 
behind that analysis is to be able to predict a-priori which of the subsystems in a given MD problem are 
“loosely coupled” in the sense that a GGS approach will be computationally efficient when applied to the 
solution of these sub-systems. The estimates allow to use mathematical theory, for example theory of Partial 
Differential Equations (PDEs) in the differential level or linear algebra in the discrete level, to estimate these 
convergence factors. Since the estimates are given in terms of the Jacobians and Hessians of the MD system 
they can be approximated by existing analysis and optimization disciplinary codes. Throughout the paper 
we assume that the problems at hand are well-posed and that there exists a unique solution; only the issue 
of numerical convergence to the solution of these problems is addressed. 

Demonstration of the convergence estimates and numerical results are given for a system composed of 
two non-linear algebraic equations and for a system composed of two PDEs modeling aeroelasticity. The 
numerical results are compatible with the analytical estimates. 

Throughout the paper the following notation is used (some of the notation is taken from [11]): 
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Ai - Analysis of system i, 

S{ - Sensitivity (Adjoint) analysis of system i, 

Oi - Optimization of system 

(Z U Z 2 ) - Single execution of Z\ followed by Z 2 , 

[Z\,Z 2 \ - Nested execution of Z\ followed by Z 2 , 

[Z 1 IIZ 2 ] - Nested execution of Z\ and Z 2 in parallel (freezing the coupling information at each execution). 

For example, {A\,A 2 ) denotes the analysis of system 1 followed by the analysis of system 2, and 
[0\(A\, S\)\\0 2 {A 2 , S 2 )\ denotes a nested execution of solving in parallel the MDO problem composed of 
disciplines 1 and 2 (freezing the coupling information at each execution), where in each disciplinary opti- 
mization problem there is a sequential application of the disciplinary analysis and sensitivity analysis. 

The paper is divided into sections as follows: 

In Secs. 2, 3, and 4 the MDA, MDS (Adjoint), and MDO problems are treated respectively. 

In Sec. 5 the theory is demonstrated on a non-linear algebraic set of two equations. 

In Sec. 6 the theory is demonstrated on a set of two PDEs modeling static aeroelastic system. 

In Sec. 7 discussion and concluding remarks are made followed by the list of references. 

Appendix A presents the norm used throughout the paper. 

Appendix B gives a short presentation of the sensitivity equation and the adjoint methods. 

Appendix C gives the proof of Lemma 1 . 

Appendix D gives the Fourier analysis of the aeroelastic problem which is discussed at Sec. 6. 

2. The Multidisciplinary Analysis Problem. In a Multidisciplinary Analysis (MDA) problem a 
coupled system of equations (typically PDEs) is solved, 

(2.1) Rj(Qi,---,Q N )=0 ; j = 

where Rj denotes the j ’th state equation and Qk denotes the k ' th state variable. If R 3 is a differential operator 
then in the discrete level each of the equations Rj is composed of a system of Nj algebraic equations and 
each of the state variables Qj is a vector of the size Nj . 

The notation (Ak x , A* 2 , • • • , Ak N ) is used to denote a GGS iteration for a MD system which is composed 
of N state equations, where the integers (fci, k 2 , ■ • • , kw) are a permutation of (1, 2, * • * , TV); equation k\ is 
solved first for the state variable Qk n then equation k 2 for Qk 2 and so on. The notation (Ak x \\Ak 2 1| • • • HA**) 
is used to denote a GJ iteration in which all the equations are solved in parallel. Note that in the GGS case 
there are N\ different sequences of iterations depending on the order of solutions while in the GJ case there 
is only one choice. 

For simplicity a two disciplines system will be considered from here on: 

/ 2 2 ) ^i(Qi,<?2)=0 

Q 2 ) = 0 . 

The system (2.2) can be solved in a tightly coupled iterative algorithm for example by Newton’s method 
or by an exact solver for the full system. One difficulty with a tightly coupled solution is that the full problem 
is typically ill-conditioned; even if each of the disciplinary problems is well-conditioned (by preconditioning 
for example), the Jacobian’s eigenvalues of one discipline will typically be different than the eigenvalues of 
the second discipline. Another practical difficulty with this approach is the need to join together the analysis 
codes of the different disciplines. 
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2.1. A Generalized Gauss-Seidel (GGS) Scheme. 

Definition: 

A GGS iteration for the two discipline system (2.2) is defined by the following procedure (see Fig.(l)): 
Scheme 1.1 

Starting with an approximation (fr° m here on the notation n will be omitted, for example 

1. Solve for Q ” +1 the equation J?i(Qi +1 , Q 2 ) = 0, 

2. Solve for Q 2 +1 the equation Q2 +1 ) = 0. 



Fig. 1. GGS Iteration for the MDA Problem: {Ai^A^). 


In some cases the GGS sequential algorithm might require many MD iterations thus it might result with 
a much higher computational cost than a tightly coupled solution. In such cases there is also a practical 
difficulty since each of the disciplinary solutions is computed with a different code, thus each GGS iteration 
requires information interchange between these codes. 

However, in some cases the GGS scheme converges efficiently to the MD solution. In the rest of this 
section the condition for convergence of GGS iteration is analyzed and a quantitative estimate of its con- 
vergence rate, fi a , is derived. This might serve as a tool to develop algorithms in complex systems where 
possibly only a subset of the system may require tightly coupled solutions; for a system composed of more 
than two disciplines it is possible to have a subsystem which is solved in a tightly coupled manner. For 
example, in a three disciplines system it is possible to define the MD iteration (Ai, [A 2 , A 3 ]), i.e., solve state 
equation number one for Q\ and then solve tightly coupled the state equations number two and three for 
Q 2 and Q 3 . It is possible that different sequences of GGS iterations result in different error reduction, e.g., 
(Ai, A 2 ,A 3 ) versus (A i? A 3 , A^. In that case we would like to analyze which of the sequences results in a 
maximum error reduction. 

2.1.1. Linearization. A linear relationship between errors and residuals for the system of state equa- 
tions is assumed. In the two disciplines system such a relationship can be written in a matrix form as 

(2.3) ( Rl ' Ql Ri ' Q2 ^ ^ = _ ^ \ 

\ R 2>Qi R 2 ,q 2 ) \ Q 2 ) \ R 2 ) 

where in general the quantities Rij depend on the state variables Q and where Q and R denote the errors 
and residuals at iteration n respectively. The block matrices /L q . are dominated by the Jacobian f#- and 
are equal to it in linear systems. We will assume that the block diagonal matrices in (2.3) are invertible. 

Note that in the differential level the quantities should be understood as the Frechet derivative of 
a function space entity, Ri € W, with respect to another, Q } e X. Such a derivative is an operator from X 
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to W: 


R 3 , Qj : X — W. 


2.1.2. The GGS Iteration in terms of the System Jacobian. Starting with the states Q\ and 
Q 2 the GGS iteration (see scheme 1.1) can be written as 


(2.4) 


Qi =Qi + Qi 
Q " +1 = Qi + Q 2 

where Q\ and Q 2 satisfy the following system of equations: 


(2.5) 



--(!)■ 


2.2. Convergence Analysis. 
Definition: 


The convergence factor for the analysis problem (2.2), /i G , is defined to be the maximum of the disciplinary 
error reduction in the state variables as a result of one application of the MD iteration: 


( 2 . 6 ) 


r ll<?? +1 l 
* 1 HOill 


IIQ2 +1 II ii^ +1 
IIQ 2 II ’ ’ 


\\Q 


o\ 

Nil y 


In non-linear problems the convergence factor depends on the state variables, fj, a = yi a {Q). 

For small values of fi a the MDA problem can be solved with a small number of GGS iterations. The 
higher the value of (but still smaller than one) the more GGS iterations required to solve the MDA 
problem. 

The convergence of the loosely coupled iteration is determined by the map T which relates the error, e, 
to the the error after application of a MD iteration, e n+l : 


(2.7) 


T : e 


e n+l . 


In the following the map T is studied for the GGS and GJ iterations for the solution of the MDA problem. 
We first find the explicit form of the map T (locally) and then take its maximal singular value, cri(T), as the 
measure of its L 2 norm (see Appendix A). The error decrease as a result of an MD iteration is estimated 


(2.8) 


ll^ +1 l 




We claim that any MD iteration, including the sequential GGS (which is characterized by a map T defined 
in (2.7)) converges to the MD solution if and only if <Ji(T n ) < (1 — e), for n > N where N is some integer 
greater than unity, and e is some real positive number smaller than unity: 0 < e < 1. 

Lemma 1; 

For a two states system (2.2), starting with errors Q\ and Q 2 in the state variables and assuming that 
the errors satisfy the quasi-linear approximation (2.3), a GGS sequential iteration (Ai,^) results in the 
following error relations: 


(2.9) 


Qi +1 = (^Q^l ,<hW&)- l KQ\)& 
QZ +l = (^Q 2 R2 , Qi R^ Q Ri,Q 3 )Q2. 


The proof is given in Appendix C. 
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2.2.1. The Convergence Matrix. The convergence matrix for a two discipline system is defined by 


( 21 °) = ^ Q R2, Qi R^ Q R^Q,- 

In terms of the convergence matrix the errors in the state variables are given by: 

(2- 11 ) Q2 + 1 = ( II ^a)02 

k=n 

( 2 i2 ) <5? +1 = -(fl ll Q 1 r 1 /i ll o,( n 

k=n — 1 

Defining a new variable for discipline 1: 

Z\ = —R 2 ,q 2 R2,QiQ\ 

we get that the error reduction in z\ is identical to the error reduction in <2 2 at all iterations and therefore 
we conclude that the convergence matrix determines the total error reduction in the system. 

The local error at iteration n is bounded by (the same relation holds for z\) 

(2-13) ||Q2ir +1 < ||A< a ||||Q a ||. 

2.3. Connection with Classical Gauss-Seidel and Jacobi Methods. The iteration (2.4-2.5) is a 
generalization of Gauss-Seidel (GS) iteration to systems of equations. When applying GS iteration to solve 
an algebraic system of equations Ax = b where A is a matrix and x and b are vectors with proper dimensions, 
then a GS iteration has the form 

(2.14) Mx n+1 = Nx + b 

where M = D + L, TV = —17, and D , L , and 17 are the diagonal, lower triangular, and upper triangular of A 
respectively. The convergence of GS iteration is determined by the spectral radius p(M~ x N) ([9] p.508): 

Theorem 1: 

Suppose 6 e JR and A = M + N e M nXn is nonsingular. If M is nonsingular and the spectral radius 
of M~ l N satisfies the inequality p(M~ l N) < 1, then the iterates x W defined by = Nx (/c) + b 

converge to x = A~ l b for any starting vector 


Generalized Gauss-Seidel 

The MD iteration defined by (2.4-2.5) is a GS iteration (2.14) with 


(2.15) 


M = 


Rl,Qi 0 

^2 ,Q X ,Q 2 


N = — 


0 R\ ,q 2 \ 

0 0 ) 


The matrix M 1 N is given by 


(2.16) M~ X N = ( ° , ^ 

\ 0 ^2,Q a R2,Qi Rl,Q 1 Rl,Q3 j 

By theorem 1 the iteration (2. 4- 2. 5) converges if p(M~ 1 N) < 1 where 
(2-17) P (m-'n) =p(^ j R 2 ,q 1 ^ iJ R 1 ,q j ). 
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Note that the expression in (2.17) is the convergence matrix (2.10). 


Generalized Jacobi 

The sequential iteration (Ai,A 2 ) can be replaced by a parallel version (^4^ ||^4 2 ) where A\ is solved in 
parallel to A 2 . This iteration will be called a “Generalized Jacobi” (GJ) MD iteration, since it is a Jacobi 
iteration with 


(2.18) 


M — 


Ri, Ql 0 \ 

0 R 2 ,q 2 j 


The matrix M 1 N is given by 


( 0 R UQ2 \ 
\ 0 / 


(2.19) 


M~ l N = 


R 2,Q 2 R 2,Qi 



According to theorem 1 the GJ iteration converges if p(M~ l N) < 1 . 

Eq. (2.19) implies that the asymptotic convergence rate of the GJ scheme is the square-root of the GGS, 
i.e., the sequential (GGS) approach requires less MD iterations than the parallel (GJ) to converge. 

In a similar manner other iterative schemes can be generalized to MD systems (c.g., Richardson, Kacz- 
marz, Etc.). Also, it is possible to apply generalized acceleration techniques developed for systems of equa- 
tions to MD systems, foe example, SOR of GS can be applied to accelerate the numerical convergence of 
the GGS scheme; in that case the update will depend on the SOR parameter, cj, and one should choose the 
SOR parameter, u, such that min^ ||Al a (u/)||. Extension of this work to acceleration of MDO schemes will 
be discussed elsewhere. 


3. The Multidisciplinary Sensitivity (Adjoint) Problem. In a Multidisciplinary Sensitivity Anal- 
ysis (MDS) problem the sensitivities of the states or of a cost functional with respect to disciplinary “design 
variables”, (61 , • • ■ ,£>m)> are computed. In the discrete level each bk is a vector of length M*. The system 
(2.1) now includes the design variables and have the following form 


(3.1) R j (Q u ---,Q N ,b k ) = 0 ; j = l,-.-,JV ; keC 3 

where Cj is a subset of the integers of { 1 , • • * , M}. 


In this section the problem of computing sensitivities of a cost functional, J, with respect to the design 
variables is addressed. The cost functional, /, is a mapping between a vector space (a Hilbert space in the 
PDE level) to the real numbers and is possibly dependent on all the variables: 

/ = I(Q !, • • • , Q n , 61 , • • • 


The “MDS problem” is the computation of the sensitivities: 

dl 


(3.2) 


(5i 


/ dl dl \ 

' ,9M ^ = VdV'"’d6‘ m)' 


As in the previous section, for the sake of simplicity, the discussion will be restricted to the following system 
of equations which is related with aeroelasticity: 


R i{Qi,Q2,h) = 0 

r 2(QuQ 2, b 2 ) — 0. 
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In the aeroelastic application Q\ stands for the flow variables, Q 2 for the structural deflection, 61 the 
aerodynamic shape and b 2 for the structural rigidity (see Sec. 5). Similar analysis can be done to any other 
MD system. 

Since there is only a single quantity as an “output” (the cost functional I) and many “inputs” (the 
design variables {&*}££ j) the most efficient method to compute the sensitivities in such a case is the adjoint 
method . A brief explanation of the adjoint method is given in appendix B (for a more rigorous treatment of 
the subject see for example [12]). 

Since the Jacobian in the adjoint equation is the adjoint of the Jacobian in the sensitivity equations 
(obtained by direct differentiation) the convergence estimates for both are identical (see Appendix B). 


3.1. Adjoint Formulation. The derivatives of the cost functional with respect to the design variables 
(sensitivity gradients) are given by 



where the Lagrange multipliers, Ai and A 2 , satisfy the following “adjoint” or “costate” equations” (If, k 
denotes ~ and R* Q . denotes (|^) ): 


Pi(Ai,A 2 ,Qi,Q 2 ,&i,& 2 ) = Jfc! + Al -f i^ Q t A 2 = 0 
P 2 (Ai, A 2 ,Qi,Q 2 ,6 i,6 2 ) = h 2 +^i,q 2 Ai + R^q 2 A 2 = 0. 

Note that the costate equations are linear in the costate variables, A, with non-constant coefficient that 
depend on the state and design variables, Q and b respectively. 


3,2. Multidisciplinary Feasible Sensitivity Solution: [Si, £2] = [Si[Aj, A 2 \, S 2 [A\, A 2 ]]. Each of 
the disciplinary sensitivity calculations is assumed to be done on a MDA solution as shown in Fig. 2. Since 
the state equations are assumed to be solved prior to the solution of the MDS problem only the system (3.5) 
is considered. In terms of the notations used in the previous section, the system of error equations is 
given by 


(3.6) 


( r Iq, 

V «Uh 


*5*1 \ ( M = 

) V ^ ) \ p * ) 


where Pi and P 2 are the residuals of the system (3.5) at iteration n. A GGS iteration and the convergence 
factor are defined in similar to the MDA case and the convergence factor is given by a similar formula. 
Notation . 

The notation for a GGS iteration for a MD adjoint problem which is composed of N costate equations 
is (Sk ! , Sk 2 , * * * , Sk N ) where the indices should be understood as in the analysis case. 


Scheme 2.1: A GGS Adjoint Iteration 

A GGS iteration of the two discipline system (3.5) is given by the following two steps, starting with an 
approximation (Ai, A 2 ) and keeping the state (QJ, Q 2 ) and the design (61, & 2 ) variables fixed: 


1. Solve for A" +1 the equation 

P 1 (A? +1 ,A 2 = A2,Q* 1 ,Qlb 1 ,b 2 ) = 0 
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2. Solve for A£ +1 the equation 

P 2 (A 1 = A^AZ + \Ql,Q* 2 ,b 1> b 2 ) = 0. 

SI S2 

a” a d 
1 2 

Fig. 2. MD Feasible Sensitivity Solution: [£i,S 2] = [5*1 [>li , A 2 ], Si[A\, A 2 ]]. 


Lemma 2.1: 

For a two states system (3.6), starting with errors Ai and A 2 in the costate variables and assuming fixed 
states, a GGS iteration of the adjoint problem (Si, S 2 ) results in the following relations for errors reduction: 

(37) A " +1 = [^Q 1 ^,0 1 (« n - 1 ) 2 “,Q a (^- 1 )t.Q 2 ]A 1 

A 2 +1 = [Kq, r Iq, r Vq^,q]^- 

The proof is identical to the proof of Lemma 1 . 

The convergence matrix for the sensitivity problem is defined by 

( 3 -8) -H,. = iq; Q ^ Ql R-^Ri Q2 . 

3.3. Single Discipline Feasible Sensitivity Solution: [*Si(Ai), £ 2 (^ 2 )]. I n this section the conver- 
gence factor for the combined computation of analysis and sensitivity, of a MD system as shown in Fig. 3, is 
estimated. 

(A.S,) (VS,) 

Qf Q? 

A. A “ 

I 2 

Fig. 3. Single Discipline Feasible Sensitivity Solution: [Si(v4i ), ^ 2 (^ 2 )] - 
Scheme 2.2: A GGS Analysis and Adjoint and Iteration 

The process which is under consideration is summarized by the following two steps iteration, starting 
with an initial approximation (Q®, A®, A°) iteration n + 1 (where n > 0) is given by: 

Step 1: (Ai,5i) , 

1. Solve for Q" +1 the state equation, 

Rx{Q n l + \Q2 = QlM) =0. 




9 





2. Solve for A" +1 the costate equation. 

Pi ( A? +1 , A 2 = A 2 ,Qi +1 ,Q 2 ,h,b 2 ) = 0. 


Step 2: (A 2 ,S 2 ) , 

1. Solve for Q 2 +1 the state equation, 
fl 2 (Qi = Q? + \g 2 + 1 , 62 ) = o. 

2. Solve for A£ +1 the costate equation, 

ft(Ai = A?, A” +1 ,£?? +1 , <?5 +1 , 6 1? 6 2 ) = 0. 


Definitions: 

The sensitivity convergence factor , /x s , is defined to be the error reduction in the costate variables as a 
result of one GGS analysis and adjoint iteration as defined by Scheme 2.2. 


Lemma 2,2: 

For a two states system (3.3), starting with errors Ai and A 2 in the costate variables and assuming 
that the errors in the state variables satisfy the quasi-linear approximation (2.3), then a GGS iteration 
((Ai, S \ ) , (^ 2 , £ 2 )) f° r the solution of the MDS problem results in the following errors reduction: 

e" +1 = S~ 1 1 S 12 (S2 2 1 )- 1 SZr 1 ei 

e" +1 = S^ 2 1 S 21 S n 1 S l2 e 2 . 

where 
(3.9) 

and the errors e* = (Qi, A*) 

Proof: 

The GGS analysis and adjoint iteration is updating the variables (Qi, A\) and (Q 2 , A 2 ) with the quan- 
tities e\ and e 2 respectively which satisfy the following equation 

(3I0) -U) - (£ £)U) 

where r» = (Ri,Pi). 

We brought the MDS problem into the same form as the MDA problem (see Eqs. (2.3) and (2.5)) with 
the changes Rij — > 5^, Rj — ► rj and Q e, and Q — > e. Using Lemma 1 we conclude that 

= ^n 1 ^ 12(^22 *) l &2i 
e 2 +1 = SgSaiSrfSuea. 

□ 


Ri 

P 


0 \ 
\Qi P'Aj ) 


10 



The convergence matrix of a GGS sensitivity (adjoint) iteration, ($i(Ai), $ 2 (A 2 )), is defined by 


(3.12) 


•Ms — $22 *$21$n 1 $12- 


The convergence factor for the iteration $ 2 (^ 2 )) is bounded from above by the maximum singular 

value of Ad 5 or equivalently by the square- root of the spectral radius of Ad* Ad 3 : 

(3.13) 


Note that the block diagonal terms in the convergence matrix, Ad 5 , are Ad a and Ad a *: 


(3.14) 


M s 


Ma 0 \ 

A4q$ A4a* J 


The eigenvalues of Ad s are equal to the eigenvalues of Ad a and Ad a *. However the maximal singular 
value of Ad s is different than that of Ad a and Ad a * and therefore the convergence rate of the scheme 
[$i(Ai), S 2 (A 2 )\ is possibly different than that of [S\[Ai, A 2 ], S 2 [A 2 , A 2 ]\. 


4. The Multidisciplinary Optimization Problem. In a Multidisciplinary Optimization (MDO) 
problem the minimum of a given cost functional with respect to the disciplinary design variables, (b\, * • ■ , 
is computed: 


(4.1) 


min ,QN,bw-,b M ), 

(61, -,OM ) 


subject to the MD state equations: 


Rj(Qu ■ * ' 1 Qn, bk) = 0 ; j = 1, ; keCj. 


As in the previous sections, for the sake of simplicity, the discussion will be restricted to the following 
problem which is related with aeroelasticity: 

(4.2a) min /(Qi, Q 2 , b u b 2 ) 

(fci , 62 ) 

subject to 


(4.2b) 


Ri(Qu Q 2 , bi) = 0 

R2(QuQ2,b 2 ) = 0. 


The adjoint formulation is used for the derivation of the sensitivity gradients (see Sec. 3.1 and Appendix 


B). 


4.1. Multidisciplinary Feasible Optimization Solution: [0i,0 2 ]- Each of the disciplinary opti- 
mization problems is solved on a MDA and MDS feasible solutions (see Fig. 4): [Oi , 0 2 ] = [0\ ([Ai , A 2 ], [$ 1 , $ 2 ]) 
, 0 2 ([Ai, A 2 ], [Si, $ 2 ])]- Since zero errors in the states (3.3) and costates (3.5) are assumed, only the system 
of the design equations (3.4) need to be considered: 

(4 3) Al + ^*2,6 i A 2 = 0 

+ R[ m Ai + ^ ib2 A 2 = 0 
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where the following notation is used h k = and R* b . = In that case the system of error equations 
for the design variables can be written as follows 



where g\ and are the residuals when solving Eq.(4.3) (these residuals are equal to the sensitivity gradients) 
and the matrix H is the transformation between errors and residuals of the design equations (gradients). 
The matrices Ha are dominated by the disciplinary Hessians (assuming feasible state solution, the first order 
term in a Taylor expansion of the gradient is the Hessian, thus, in the first order approximation, the Hessian 
relates the errors in each of the disciplinary design variables with the negative gradient of that discipline). 

Assuming MDA and MDS feasibility at each optimization step, the GGS iteration is defined as in the 
MDA case with O replacing R and where equation (4.4) replaces (2.3) (the Hessian blocks, H i3 , replace the 
Jacobian blocks, Rij). 



01 02 


min I (b , b = b“) 

0+1 

min I (h = fi'Vb ) 

b, 1 2 2 

C.,b,A 

b 2 

W<V=° 

Q,b D ,A 2 


^(A,A)=0 

^2 i 2 

FyA,A 2 ) = o 




i 2 2 


Fig. 4. MD Feasible Optimization Solution: [0i,02]- 


Lemma 3.1: 

For a two disciplines optimization problem (4.2a-4.2b), starting with errors bi and b 2 hi the design 
variables the MD feasible scheme, [Oi,0 2 ]i results in the following expressions for errors reduction: 

(4 6” +1 = [Hn 1 H 21 (H^ 1 )^H- 2 l ]b 2 

K +l = [H 22 1 H l 2 H- l 1 H 21 ]b 2 . 

The proof is similar to the proof of Lemma 1. 

The convergence matrix in that case is given by 

(4-6) M h = 

and the convergence factor is bounded by 

(4.7) fi H < p 2 

4.2. Single Discipline Feasible Optimization Solution: [0i(Ai, Si), 0 2 (A 2 , Sa)]. I* this section 
we estimate the convergence factor for the solution in which there is decoupling of the analysis, the sensitivity, 
and the optimization solutions of one discipline from the other as shown in Fig. 5. In terms of the notation we 
have used so far the decoupled scheme is denoted by [O kl (A kl , S kl )> O k2 (A k2) S k2 )] y i.e., optimizer of system 
is using repetitive calls to the analysis and sensitivity analyzers of system A;i (and uses fixed values of 
states and sensitivities of system fc 2 ) followed by optimizer of system which acts similarly for system A; 2 . 
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O (A , S ) O (A , S ) 

iii 2 v 2 r 


min I (b , b = t> ) 

□+1 n+1 n+1 

min I (fy - 6 + .'b ) 

b, ' 2 2 

^.b.A 

b2 ' 2 

*i(Q,.Q = <|) = o 

ODD 

R(Q = Q.q> = 0 

P^,A f /f) = 0 


p 1 (A 1 =A n : , A) = o 


Fig. 5 . Single Discipline Feasible Optimization Solution: [Oi (v 4 i , Si), 02 (^ 2 , £ 2 )] . 

The MD scheme under study is summarized by the following two steps scheme. Starting with an initial 
approximation, (Q? , Q®, A ?, A§, iteration n + 1 (where n > 0) is given by: 

Scheme 3.2: A GGS MDO Iteration: (Oi{A \ , S\), 02(^2, 52)) 


Step 1: Oi[j4i,5i] , 

Solve for 6" +1 the optimization problem 

min fc „ + i I(Q? +1 ,Q 2 = Q^,b^ + \b 2 = b%) 
subject to 

^i(o? +1 ,Q 2 = g?,6? +1 ) = o. 


Step 2: 02^2) 

Solve for 6£ +1 ths optimization problem 

min ^/(Qi =Q? + 1 ,gj + 1 , 6 1 = 6 ? + \& 2 +1 ) 

subject to 

R2(Qi = Q? +1 ,QP\I% +1 ) = 0 . 


Lemma 3.2: 

For a two disciplines optimization problem (4.2a-4.2b), starting with errors b\ and 62 in the design 
variables and assuming that the errors in the states satisfy the quasi- linear approximation (2.3), then a 
single discipline feasible GGS iteration, (Oi(Ai, Si), 02 (^2 > ^2))j results in the following errors reduction: 


e ” +1 = (On'OMO^'r'Oli 1 )^ 

e 2 +1 = (O 22 1 0 21 0r 1 1 0 12 )e 2 

where the matrices 0^ are defined by 



( Ri.Q, 

0 

Ri,bj 

(4.8) 

Oij = Pi,Q, 

A, 



V 9i,Q 3 

9i,A.j 

9i,bj 


and ei = (Qi,Ai,bi). 


Proof: 
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In terms of the matrices Oij the relation between the errors and the residuals is given by 

( 4 . 9 ) - ( ri ) = ( °" 0,2 V V 

\ r 2 J \ O21 O22 ) \ e 2 J 

where r t = (/?,, P,,g t ). 

The GGS MDO iteration is updating the variables (Qi,Ai,bi) and (Q 2 , A 2 , b 2 ) with the quantities e\ 
and e 2 respectively which satisfy the following equation: 

( 4 . 10 ) - ( r ' ) = ( ° ) ( !' ) . 

\ r 2 J \ O21 O22 J \ 62 J 

We brought the MDO problem into the same form as the MDA problem (see Eqs. (2.3) and (2.5)) with the 
changes — ► Oij , Rj — ► Q — > e, and Q e. Using Lemma 1 we conclude that 

(4 n) ei +1 = (ori'OntOJr 1 )-^- 1 )^ 

er 1 = (o 22 l 0 21 O u 1 0 12 )e 2 . 

□ 

The convergence matrix in this case is given by 

( 412 ) Mo = 0 22 1 0 21 0 n 1 0 12 . 

As a result of one iteration ( 0 \ (A\ , Si). 0 2 (/l 2 . <S 2 )) the total error in the variables of discipline 2 (and 
also in the transformed variables of discipline 1 (see Sec.2.2.1)) e— (Q y X,b) is reduced by 

( 4 . 13 ) e” +1 = M 0 e. 


An upper bound for the convergence factor, /z 0) is estimated by the following formula: 

( 4 * 14 ) Mo < 

5. Test Case I: a Non-Linear Algebraic Example. In this section we treat a set of two non-linear 
algebraic equations. The aim here is to demonstrate how the MD convergence estimates are computed when 
the state equations are non-linear. In this example the estimates will be computed in the same space as the 
numerical computation is performed, in contrast to the next section in which the estimates are approximated 
in the infinite dimensional space while the equations are solved in finite dimension. In non-linear problems 
the equations are first linearized and the convergence theory is applied to approximate the MD convergence 
rate on the linearized system. 

5 . 1 . Multidisciplinary Analysis. Let us examine a simple set of non-linear algebraic equations: 

^5 x \ X 2 = b\ 

= b 2 

where we set b\ = 62 = 2. A solution of the above system is given by (ari = 1,X2 = 1). 

A GGS iteration for the iterative solution of the system ( 5 . 1 ), (Ai, A 2 ), is the Gauss-Seidel (GS) iteration. 
The convergence matrix (a scalar in this case) for the “MD” iteration (Ai, A 2 ) is given by (Eq. ( 2 . 13 )): 

< 5 - 2 ) Ma(« > 1) < ^Q 3 R2, Ql Rl l Q Rl,Q 3 = T-Lg. 


14 



In this case the convergence factor depends on the states, X\ and x 2 . At the solution, (xi = l,x 2 = 1), the 
convergence factor is asymptotically equal to | thus the predicted asymptotic convergence rate for the GS 
iteration for the solution of Eqs.(5.1) is fi a ~ 

Apphcation of formula (C.8) for the first iteration results in the following inequality 


(5.3) 


M a(n = 


1) < max 


r_l_Jx|i 1 1 

1 3xf |x°| ’ 9x^X2 J 


Fig. 6 depicts the actual convergence rate versus the theoretical upper bound, where ( 1 theory given in 
(5. 2-5. 3), and 


factual — max 


r l*i +1 l l^2 +1 h 

i 1*1 I ’ 1^2 1 i 


For n >2 both /i theory and factual are equal to //, = 


5.2. Multidisciplinary Adjoint. We define the cost functional 

(5.4) I(xi , x 2 ) = (xix 2 - l) 2 + 2 X 1 + 

where x\ and X 2 satisfy Eq.(5.1). 

The sensitivities ^ and ^ are given by 


(5.5) 


9 1 


dl _ 
dbi ~ 
^ _ d/ _ 

22 - db~ 2 ~ 


— Ai 

-a 2 


where the costates, Ai and A 2 , satisfy the adjoint equations: 


(5.6) 


3xf Ai + A 2 = -I x 1 = “ ^2x 2 (xix 2 - 1) + 

Ai + 3x^A 2 = -I X2 — - ( 2 x 1 (xix 2 - 1) + 2x 2 ^ . 


5.2.1. Multidisciplinary Feasible Sensitivity Solution: [Si[Ai, A 2 ], S 2 [Ai, A 2 ]]. Since all the op- 
erators in this case Ri,Q } are real scalars the convergence estimates are identical to the MDA case. In the 
numerical test the MDA problem is solved first and then the GS scheme is applied to the adjoint system 
(5.6). Fig. 7a depicts the actual versus the upper bound convergence rate for this scheme. As in the MDA 
case there is a good agreement between the theoretical upper bound and the actual convergence rates. 


5.2.2. Single Discipline Feasible Sensitivity Solution: [Si(Ai), S 2 (A 2 )]. We apply Lemma 2.2 to 
estimate the upper bound for the convergence rate. In this case the matrices S{j are given by 


(5.7a) 


(5.7b) 


where 


(5.8) 


S n 


3x 


2 

1 


6xiAi -F 




> 

) 


— 2X2X2 d - 2 

I Xl , X3 = - 2 

4xix 2 2 

/x2,x2 — 2xiXi -f- 2. 
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The upper bound for the convergence rate is estimated by 

Ms ^ 11^22 ^21 ‘S 11 *^12 1 1 * 

In the numerical test scheme 2.2 was applied [5i(Ai) ? £ 2(^.1 )]• Fig.7b depicts the actual versus theoretical 
upper bound of the convergence for this case. The theoretical upper bound for the first GGS iteration is 
Vtheory(n = 1) ~ 2.6 (n = 1 in the figure). The actual convergence rate approaches the MDA convergence 
rate fi a = ~ for n » 1. 

5.3. Multidisciplinary Optimization. The optimization problem is defined by 


(5.9) 


min I(x 1 , 3 : 2 ) 


subject to Eq.(5.1), where the cost functional I{x 1 ,^ 2 ) is defined in (5.4). 

5.3.1. Multidisciplinary Feasible Optimization Solution: [Oi,C?2]* The Hessian is computed 
with a combination of the adjoint and direct- differentiation method (the combination of the adjoint and 
sensitivity methods to compute the Hessian was originally introduced by R. Haftka [13]). The Hessian is 
given by 


(5.10) 


where A = |^ satisfies 


H= dg dg dQ | dg . 
db dQ db dX 


dP dPdQ dP . 
db + dQ db + d\ X 


_ / dR \ * t / dl \ * 

= ( dQ ) A+ (ag) 


(5.11) ^- + ~^ + ^-A'=0 
and P = P(Q, A, b) denotes the residuals of the costate equation 

(5.12) 

For the problem (5.9) 

(5.13) 

Thus the Hessian is given by 


— &R- — — 0 and ^ 

db dQ db u 4111(1 dx 


(::)■ 


/ H n 

h 2 1 \ 

= a' = -(%?) 

\ h 12 

h 22 / 

[dxj 


-1 


dP dQ _ 
dQ db ~ 


(5.14) 


3xf 1 

1 3^2 


-1 


6xj Al + I Xli Xi Ix\yX2 

6X2X2 ~f“ I x 


\(Zx\ 1 \ 

/ \ 1 3X2 ) 


-1 


where I Xi>Xj are given in (5.8) 

The convergence factor jin for this scheme is estimated by 


(5.15) 

where Mh is given by 

(5.16) 


Mh < P 2 ( M* h Mh ) 


Mh = 
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In the numerical test each of the sub-optimization problems was solved at a time, i.e., 


0 \ = min J(xi) ; 0 2 = min J(x 2 ), 

b i 62 

where full MDA and MDS feasibility were maintained at each solution. Fig. 8a depicts the actual versus 
theoretical upper bound of the convergence for this case. The actual and theoretical upper bound convergence 
factors for the first GGS iteration are given by /i ac * ua z( n = 1) ~ 2.7 and Htheor y {n = 1) ~ 9.5 respectively. 
The asymptotic convergence rates are given by /i ac tua/(™ = 1) « 0.93 and fitheory(n = 1) ~ 0.95. 

5.4. Single Discipline Feasible Optimization Solution: [Oi{Ai, Si), 0 2 (A 2 , £ 2 )]* The conver- 

gence of the GGS scheme, [Oi(i4i, S\), 0 2 (A 2 , S 2 )], is determined by the matrices O tJ defined in Eq.(4.8). 
For the problem (5.9) the matrices Oij are given in the following: 




( R UQi 

0 

■Rl,6i 

\ 



3xf 

0 


(5.17a) 

O 

II 

P l,Qi 

PiM 

PiM 


— 


6a?iAi ~h I x i,ii 

3x? 0 



K 9i,Qi 

9 1,A, 

9iM 

J 



0 

-1 

0 J 



R i ,Q 2 

0 

R l,b 2 \ , 

( 1 

0 

°\ 


P l ,Qi 

PiM 

PiM = | 


1 

0 


9\,Qi 

Q\^2 

9i, b 2 J ' 

^ o 

0 

0 ) 




( 

R 2,Qi 

0 

R 2M \ 

( 1 0 

°\ 

(5.17c) 

02 1 = 


P 2 ,Ql 

P2,X 1 

P2M = 

1 ^12, Xl 1 

0 



\ 

92.Q 1 

92M 

92 , bl J 

V 0 0 

0 / 




f R 2,Q 2 

0 

P 2 ,b 2 

\ 

/ 3xl 

0 

-i\ 

(5.17d) 

o 22 — 

P 2,Q 2 

P 2 M 

P 2 M 


= 6X2^2 + ^X2,X 2 

3x 2 0 



\ 92, q 2 

92 , \ 2 

92, b 2 

) 

V 0 

-1 

0 J 


An upper bound for the convergence factor is estimated by the spectral radius of A4qMo, where 
A4o = 0 2 - 2 1 0 2 i0r 1 1 0i 2 : 


Ho < pHMqMo)- 

Fig. 8b depicts the actual versus theoretical upper bound of the convergence for this case. The actual and 
theoretical upper bound convergence factors for the first GGS iteration are given by Hactuai(n = 1) ~ 1.9 and 
Htheory{n = 1) ~ 4.3 respectively. In the numerical test the solution of the problem was obtained after the 
third iteration in which the error was reduced by a factor of approximately 1.75 ■ 10~ 6 (the optimal solution 
is given by x\ — l,x 2 — = l|, and 6 2 = l|). It is interesting to note that in this case the solution 

required much less optimization steps than in the previous case where full feasibility was maintained. This 
was predicted by the analysis although the upper bound for the convergence is higher than the one found in 
practice: for n > 1, Htheory ~ 0*25. 
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6. Test Case II: A Coupled System of PDEs Modeling Static Aeroelasticity. In this section 
we treat a set of two coupled PDEs with the aim to demonstrate how estimates in the PDE level can be used 
to approximate the different convergence factors in the problem. Since the actual numerical computation is 
done in a finite dimensional space, the continuum estimates holds only in the limit of the discrete mesh-size, 
approaching zero. Nevertheless, the continuum estimates give a good approximation of the actual discrete 
convergence (as long as the discrete mesh is not too coarse). The coupling estimates can be derived also in 
the discrete level where the different Jacobians are large matrices (rather then operators). 

In this example we use Fourier analysis and Parseval’s relation to estimate the norms of the convergence 
operators involved, however it is not necessarily the only way to get such estimates. The problem under 
consideration is linear in the state variables with non-constant coefficients. However, the optimization 
problem is non-linear since the design variables multiply the states. In the general case where the state 
equations are non-linear the first step in the analysis is to linearize them resulting in a linear equation with 
non-constant coefficients similar to the problem under consideration in this section. 

This problem contains some of the physical features of static aeroelasticity. Aeroelastic optimization 
is a major MDO problem in aerospace industry (e.g., [14] and references their-in). The aerodynamic state 
equation models the flow around a body and the design variables are some parameters defining the shape 
of that body. The aerodynamic cost functional is typically the drag or some measure of the closeness to 
some specified pressure distribution. The structure state equation typically results from a finite element 
analysis of the structure and the design variables are some parameters which determine the stiffness of the 
structure. The structure cost functional is typically the weight and there are constraints on the stress on 
different elements in the structure. There is not a unique way to define an MDO problem in that field and 
hence we will assume a simple model problem that contains some of the physics of the real problem. 

6.1. Problem Definition. We consider a two dimensional potential flow over a one dimensional beam. 
The MDO problem is to compute the set of rigidity parameters, and the set of shape parameters, 

{ a ?}£Li> such that the following cost functional is minimized (we denote the discrete quantities by a super- 
script h), 

N 2 

(6.1) F h (<f> h , W h , a h , D h ) = 7 i far# - /*'*) + 72 (itf)* + 73 (Sr#)W? + 745(0?)* 

*=1 1 

subject to the inequality constraint 


(6.2) Dl> > D min 1 < i < N 

and to the following finite difference equations, 

(1 - Ml)d* x <t>* k + =0 1 < i,k < N 

a*# 1 = #(# + w t h ) 1 < i < N 

(6.3) = #* 1 < k < N 

•AjV.fc = &N- l,k 1 < k < N 

=0 1 < i < N 


(6.4) 


S£*(£?3*W) = ~Pi ((#)<) 1 <i<N 

wf = Wfc = 0 

d^Wt = dxx^N = 0 for a simply supported beam 
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where d xi d z ,d xx and d zz denote finite difference operators for the first and second derivatives in the x and 
z directions respectively. 

The pressure, p, depends on the potential by the Bernoulli law (we assume d x <f) U 0 c ): 

( 6 - 5 ) P - Poo + pooU^dtj). 

In the numerical test the far field parameters were set to unity except the far field Mach number which was 
set to zero (the latter was chosen to simplify the analytical estimate of the coupling): 

Poo poo ~ 1 } Moo = 0. 

The forth term in the cost functional was added for uniqueness of the optimal solution a* (since only the 
derivative d x a appears in the equations). The weights in the cost functional were determined to establish a 
significant coupling between the disciplines, i.e., at the optimal solution the deflection, W * , and the shape, 
a * , have the same order of magnitude: 

7i = l ; 72 = 0.01 ; 73 = -1 ; 74=1- 

In order to avoid singularities in the beam equation, D h = 0, an inequality constraint has been applied 
on the rigidity: 

(6.6) D(x) > D rniri (x) 0 < x < 1 

where the minimal value of the rigidity was set to D m i n (x) = 10 -4 . 

The computational grid consisted of an (N x N) grid on which the potential equation was solved in the 
whole domain while the beam equation was solved on the boundary ( z — 0). On each grid point on the 
boundary, 1 < i < N, two design variables were defined: and D 

6,2. MDA. The MDA problem is to solve the system of discretized PDEs (6. 3-6. 4) for <f> h and W k 
where the design parameters, a h and D h , are fixed. A GGS iteration is composed of two steps; solving 
Eq.(6.3) for <f> h keeping W h fixed, followed by the solution of Eq.(6.4) for W h keeping <j) h fixed. By Eq.(2.13) 
the convergence of the above GGS iteration is determined by the norm of the convergence matrix , M a = 
^Q2 R ^Q^ R i,Qi Rl ^' In the following Fourier analysis in the continuum level is used to estimate the 
norm of M a . Since the estimate is derived in the continuum level and the numerical test is done in the 
discrete level, the derived estimate is only an approximation (that is expected to be more accurate with the 
refinement of the discretization). 


An Estimate of the convergence factor using ParsevaPs Relation 

We now make use of Parseval’s relation to estimate the norm of the convergence matrix, A4 a , in the 
PDE level. Parseval’s relation states that the L 2 norm of the error in real space is equal to the L 2 norm of 
the error in Fourier space: 

(6-7) |j£|U» = II4IU- 

By Eq.(2.13) 

(6.8) e? +1 = M a e 2 
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and therefore (see Appendix A) 


(6-9) ||e^ +1 ||<ma X (^;>( 0 )^||e 2 ||. 

An upper bound for the convergence factor, fa , is estimated with 

(6.10) fi a < max (m^Mo) 2 = max ( - — ftkir) \ 2 

^ k \ a ) k \G* W (D)G W (D )) 

6.2.1. Numerical Results. In the numerical test the scheme [A i, A 2 \ was apphed. A random number 
generator was used to give the initial values to the state variables, <f) h and W h . At each iteration the single 
discipline analysis problem (the the potential equation for <j> h ) was solved to machine zero accuracy, keeping 
the other discipline state variables ( W h ) fixed. This was followed by solving the second discipline analysis 
problem (Beam equation for W h ) while keeping the updated values of the first discipline (<f> h ) fixed. The 
equations were solved with a second order accurate finite difference iterative method. 

Fig. 9 depicts the actual convergence rate versus the theoretical upper bound, where /x theory given in 
(6.10) and 


factual — max 


f ll<T +1 ll ll^ n+1 h 
\ ||^»|| ’ ||^"|| /■ 


Since Gw{D) in (6.10) depends on the values of the rigidity, , the estimate (6.10) was computed by 

the analysis code for the actual values of the rigidity at each iteration. An alternative could have been a good 
a-priori estimate of the minimum value of Gw(D)- The results show good agreement with the predictions 
of the convergence estimates. It is interesting to note that although the estimate for the convergence factor 
depends on the non-constant coefficient D h still the predicted values of fa are fixed (for n > 1) on the value 
fi theory ~ 7 x 10 -3 , possibly due to the determination of the extremal values of D h , d x D h and d xx D h during 
the early stage of the computation. The actual convergence rate is about factual « 1 — 3 x 10~ 3 , (lower 
than the upper bound estimate). After the fourth iteration the errors in the state variables reach the level 
of machine accuracy and the computation is terminated. 


6.3. MDS. In the MDS problem the sensitivities on the cost functional, (6.1), are computed with 
respect to the design variables D h and a h . We choose to use the adjoint method (see appendix B) to 
compute these sensitivities since it is the most efficient method if there are less “outputs” (in this case the 
only output is the cost functional) than “inputs” (2iV design variables) in the problem. 

The costate equations for that problem, in the continuum level, are given by (the costate flow and 
structure variable are denoted by Ai = A and A2 = 77 respectively): 


Flow Costate Equation 


( 6 . 11 ) 


AA — 0 

= d x r] + 271 (d xx <f> - d x f*) + 7 3 d I W' 
d x X = 0 
d x X = 0 
A = 0 


in = {0 < x < 1 ; 0 < t/ < 1} 
on z ~ 0 
on x = 0 
on x = 1 
on z = 1. 
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Structure Costate Equation 


( 6 . 12 ) 


d xx (Dd xx r]) = —d x X - 73*9x0 on 2 = 0;0 <x<l 

* 7 ( 0 ) = * 7 ( 1 ) = 0 

d xx r](0) ~ dxxV(l) = 0 for a simply supported beam . 


Having defined the states and costates, the sensitivity gradients, g\ = and 92 = are given by 


(6.13) 


9 i = d x \ + 74 a 
92 = d xx rfd xx W + 


on z = 0 ; 0 < x < 1 . 


The MDS problem, for this example, is to solve Eqs.(6.11-6.12) for A and 77. We apply two schemes: 
the MD feasible sensitivity solution [Si, £2] and the single discipline sensitivity solution [Si(Ai), £2(^2)] 
requiring only single discipline state feasibility at each step (see Secs. 3. 2-3. 3). 

Multidisciplinary Feasible Sensitivity Solution: [5 1? £2] 

By Lemma 2.1 (Eq.(3.8)) the convergence matrix for the MD feasible sensitivity scheme is given by 


(6.14) 

We conclude that 

(6.15) 


M a * — ^2,Q 2 ^2,Qi^1,Q 1 ^1,Q 2 “ 


2kn 


G* W {D) 


Ml 


\\M a .\\ = \\M a \\, 

thus we predict the same convergence factor as in the MDA case. 

Single Discipline Feasible Sensitivity Solution: [Si (Ai), £2(^2)] 

The symbol of the convergence matrix for that scheme is given by (see Eq.(3.14)) 


•M-s — £22 £2 1*^11 £12 


-it 


M a 


0 


where 

(6.16) 


M 


«= (m,*)* = 


2k7Y 


4/C7t( 73 — 7i2/:7r) 


Mas = 


Gw(D) 


G W {D) 

An upper bound for the convergence rate is estimated by 

(6.17) ( jl s < maxpi ^A^*(fc)A4 s (A;) j . 

6.3.1. Numerical Results. In the numerical test the two schemes [£1 , £2] and [£i(Ai), £2(^2)] were 
tested. 


Multidisciphnary Feasible Sensitivity Solution: [£1, £2] 


In the numerical test the MDA problem was solved first and then the GGS scheme was applied to the 
adjoint system (6.11-6.12). A random number generator was used to give the initial values to the costate 
variables: X h and rf 1 . The same procedure as in the MDA case was repeated for the adjoint equations. 
The theoretical upper bound in this case is identical to the one in the MDA problem (Eq. 6.15). Fig.lOa 
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depicts the actual versus the upper bound convergence rate for this scheme. For iteration number n > 1 
the theoretical upper bound settles on /z theory ~ 8 x 10“ 3 while the actual convergence rate for n > 1 is 
factual ~ 2 4 X 10 


Single Discipline Feasible Sensitivity Solution: [Si(Ai), 52(^2)] 


In this case a random number generator was used to give the initial values to both the state and costate 
variables: <j > h , W h , \ h and t] h . At each iteration the single discipline analysis and adjoint problems were 
solved sequentially (the analysis problem is independent of the adjoint problem). 

Fig. 10b depicts the actual convergence rate versus the theoretical upper bound, where /z theory given in 
(6.17) and 


factual 


= max 


{ 


]|0n+i + A w ~^i || 
\\4> n 4- A n || 5 


||PF"+i+77 n+1 |h 

||^ n + fp|| /’ 


For iteration number n > 1 the theoretical upper bound, computed by (6.17), settles on /z theory ^l.lxlO -2 
while the actual convergence rate for n > 1 is factual « 1 — 5 x 10 -3 . 


6.4. MDO. The MDO problem is to solve the optimization problem (6. 1-6.4). The MDO problem is 
divided into two sub-optimization problems: 0\ and 02* In 0\ the design variable is a h while D h is fixed 
and in O 2 it is visa versa. In the numerical test two schemes were applied: the MD feasible scheme, [Oj, 0 2 ], 
and the single discipline feasible scheme, [0\{Ai, Si), 02(^2, S2)]. 

6.4.1. Fourier Analysis for the MD Feasible Scheme [Oj, O2]. By Eq.(4.6) the convergence matrix 
in this case is given by 


(6.18) 


M H = ff 22 +21+,+12. 


Following [10] the symbol of the Hessian for this problem is computed by the following procedure. The map, 
in Fourier space, between the errors in the design variables and the state variables is computed: 


/ <j>{k) \ _ / 27 rk —2nik \ / — 27riA: 0 \ / d(Ar) \ 

V W(k) ) ~ \ 2nik G W (D) ) \ 0 G d (W) ) \ D(k) ) 


The solution of (6.19a) is substituted in the following equation for the map, in Fourier space, between the 
errors in the design variables and the costate variables: 

( A(A:) \ _ / 2nk -2nik \ / 27i(27rA;) 2 -73(2^) \ / <j>(k) \ 

V *K fe ) ) ~~ V 27 Tik G W (D) ) \ 73(2tt ik) 0 ) ^ W(k) ) 

(6.19b) -( m , 2 ’ ik y '(° . ° w v 

V 27rik ) V 0 Gd fa) ) \ ) 

Finally, both the solutions of (6.19a) and of (6.19b) are substituted in the expression for the symbols of 
the sensitivity gradients to compute the map, in Fourier space, between the errors in the design variables 
and the residuals of the design equations (the gradients): 
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(6.19c) 


f 9i(k) \ _ / 0 0 \ / 4>(k) \ f 2 t Hk 0 \ / X(k) \ 

V g 2 (k) ) \ 0 G d ( v ) ) ( W(k) ) \ 0 G d (W ) ) ( f)(k) ) 

/ 74 0 \ / a ( k ) \ 

\ 0 -72 \D-l ) V D(k) ) ■ 

As discussed in [10] the symbol of that map is the symbol of the MD Hessian: 

(6 20) ( ° l{k) ) = -( Au{k) ^ 12W ) ( &{k) ) 

V 92 (k) ) V An (k) H 22 (k) ) \ D(k) ) ' 

Having computed the symbol of the convergence matrix for that problem, 

M H (k) = H^(k)H 21 (k)H^(k)H 12 (k), 
an upper bound for the convergence rate is estimated by 
(6.21) g, H < maxp^ ^M* H (k)M H {k)j . 


6.4.2. Fourier Analysis for the SD Feasible Scheme \Oi(Ai, Si),0 2 (A 2 , S 2 )]. By Lemma 3.2 the 
convergence matrix in this case is given by 

(6-22) M 0 = 0^0 21 0^0 12 

where the operators Oij are defined in (4.8). 

The matrices O tJ are given in the following: 



/ &i,Qi 

0 

Rim ^ 

{ ( 2nk 

0 

—2nik \ 

(6.23a) 

C» 

II 

>> 

b 

A, A, 

PiM 

I = [ 27 i( 2 ?rfc ) 2 

27rAr 

° 


V 9hQi 

9i,\i 

9iM ) 

/ V o 

2nik 

74 / 




( A,0a 

0 

A ,b 2 ^ 




—27 rik 

0 

0 

( 6 . 23 b) 

0\2 = 


A ,q 2 

A,a 2 

A, b 2 


= 


-73(27 rik) 

— 27 riA: 

0 



V 9i,Q2 

9i ■> A 2 

5l,6a ) 




0 

0 

0 




f A,Qi 

0 

^2,bi 

) 


1 

f 2irik 

0 

°\ 

( 6 . 23 c) 

021 ~ 


A,Qi 

A, A! 

A, 61 

1 

= 

' 

7 3 ( 27 rifc) 

27rik 

0 




V <72,0! 

9 2 , A! 

92, bi 

) 



< 0 

0 

0 / 




( A,Qa 

0 

# 2,62 

\ 

( Gw(D) 

0 

G d (W) \ 

(6.23d) 

022 ~ 

1 A ,02 

p2,A 2 

A, 62 



Ow{D) 

Gd{i 7 ) 



\ 92,02 

$2, a 2 

^ 2,62 

) 

V G d (tj) 

G d {W) 

-72 \D-i ) 


Using Parseval’s relation an upper bound for the convergence factor is estimated by the maximal singular 
value of Mo, or equivalently by 

(6.24) < ma xp* (M* 0 {k)M 0 (k)), 

k 

where M 0 — O^^iO^O^. 


23 



6.4.3. Numerical Results. In the numerical test the two schemes [0 \ , 0 2 ] and [0i(Ai, S x ), 0 2 (A 2 , S 2 )] 
were tested. 


MultidiscipUnary Feasible Optimization Solution: [0i,0 2 ] 


In the numerical test the MDA and MDS problems were solved to machine zero accuracy at each 
optimization step. The theoretical upper bound for the convergence rate was computed by (6.21), and the 
actual rate by 


factual 


rjlg^ll ll£ n+1 lh 
1 ||a»|| ’ ||D»|| I 


In that case the error in a k reduced at the first iteration by a factor of 1.3 x 10~ 5 while that of D h by 3.7 x 10 -6 . 
Since the error reduction is so significant the computation is terminated after two iterations due to the low 
numbers involved in the computation at this stage. The theoretical upper bound is fi theory ~ 3.5 x 10 ~ 2 
which is much higher than the actual convergence rate. 


Single Discipline Feasible Optimization Solution: [0i(Ai, Si), 0 2 (A 2 , S 2 )] 


In this case at each iteration the single discipline optimization problem was solved while keeping the 
variables of the other discipline fixed. 

Fig. 11 depicts the actual convergence rate versus the theoretical upper bound, where jj, theory given in 
(6.24) and 

_ f Il0 n+1 + A n+1 4- a n+1 || \\W n + l +fj n + 1 + D n+1 1| 1 

factual - max | ||^n + A" + d n || ’ || W n + fj n + D n || /' 

In this case the actual convergence rate for iterations 2,3,4, and 5 are given by 5.5 x 10 -4 , 1.84 x 10~ 3 , 
8.8 x 10” 4 , and 1.02 x 10 -2 . The theoretical upper bound is iitheory ~ 3.2 x 10 -2 . As discussed in [10] the 
fast convergence of sequential algorithms applied to the aeroelastic model problem indicates of weak coupling 
between the two disciplines. 

7. Discussion and Concluding Remarks. We present a quantitative analysis to determine the 
convergence rate of loosely coupled schemes for solving MD systems of equations. 

The sequential scheme, where the problems are solved in the disciplinary level “one- at- a- time” , is iden- 
tified as a generalization of the Gauss-Seidel (GGS) iterative method. The parallel version of the sequential 
scheme, where the problems are solved in the disciplinary level in parallel exchanging coupling data after 
each MD iteration, is identified as a generalization of the Jacobi (GJ) iterative method. Other iterative 
schemes can be generalized to MD systems (e.g., Richardson, Kaczmarz, Etc.). 

The asymptotic convergence rate for the GGS scheme can serve a measure for the coupling between the 
disciplinary problems since the looser the coupling is the lower is the value of that convergence rate. An 
upper bound for the convergence rate (convergence factor) is estimated for the multidisciplinary analysis 
(MDA), sensitivity (MDS), and optimization (MDO) problems. 

Five different basic schemes are considered (given here for problems composed of two disciplinary sub- 
problems). In all the schemes the convergence matrix is given in the general form M = Z 2 iZ 1 “ 1 1 Z X2 . The 
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quantities Z^ depend on the scheme under consideration and are summarized here for the different schemes 
presented in the paper: 


• Sequential GGS scheme for the MDA problem, [Ai, A 2 ]: Zij = Ri,Q 3 

• Sequential GGS scheme for the MDS problem 

- MD feasible sensitivity solution, [5 i,S 2 ]: Zi sJ = R* q. 

- Single discipline feasible sensitivity solution, [5i(Ai), 5 2 (A 2 )]: 

2 „ - ( ° ) 

V P,,Q, P„ A, J 

• Sequential GGS scheme for the MDO problem 

- MD feasible optimization solution, [Oi, 0 2 ]: Zi j — Hij 

- Single discipline feasible optimization solution, [0i(Ai, Si), 0 2 (A 2 , S 2 )]: 

/ 0 Ri,bj \ 

Zi,j = I Pi ,Qj Pi,Aj Pi t b 3 I ■ 

V 9i } Qj 9i,Aj 9i,bj / 


For all of these schemes an upper bound for the convergence factor is estimated by ji < <Ji(.M), where a\ 
denotes the maximal singular value of M. In problems governed by PDEs the quantities M are operators, 
typically non-differential, which measure the coupling between the disciplines in the problem. 

The convergence theory was applied to two test cases: a system of two non-linear algebraic equations 
and a system of two PDEs which model an aeroelastic system (in both test cases there are non-linearities 
causing the coupling to be a local property). In both test cases MDA, MDS, and MDO problems were 
defined, analyzed, and tested numerically. In the first test case it is demonstrated how the convergence 
factors can be estimated in a finite dimensional level. In the second test case it is demonstrated how the 
convergence factors can be estimated in the PDE level; such estimates are expected to be more accurate 
as the discretization refines. The symbol of the convergence operators were computed explicitly and where 
estimated on the boundary (which is the interface between the two disciplines involved). The symbols 
indicate that the convergence operators are pseudo- differential and contain the information of the coupling 
between the two disciplines. We think that these operators play a major role in the theory of MDO governed 
by PDEs. The numerical results are compatible with the analytical estimates. 

In applications, where large codes are involved, the different convergence factors can be approximated 
numerically, locally, using numerical computation (in the disciplinary level) of the Jacobians and Hessians of 
the different disciplines in the problem. These convergence factors can serve as a tool to develop algorithms 
for the solution of MD problems. When a problem is composed of more than two disciplines the convergence 
factor can be estimated between all the different disciplines. For example, if the problem is composed of 
three disciplines and the convergence factor of the iteration (Ai, A 2 ) and that of (Ai, A3) is much smaller 
than that of (A 2 , A3), then the scheme (Ai, [A 2 , A3]) is expected to have desired convergence properties. 

For MDO problems which have a “system level”, e.g., for 0[(Ai, «Si),(A 2 , 5 2 ),(A 3 , 5 3 )] the convergence 
theory can be applied to define algorithms for the solution of the MDA and MDS problems. Extension of 
this work to acceleration of MDO schemes will be discussed elsewhere. 
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Appendix A. Norms. Throughout the paper the L 2 norm is used. In the following we define it for a 
vector, a matrix, a function, and a differential operator: 

• Let u be a vector in C n , then 


Let A be a matrix in C mXn , then, 


||A|| = sup 


\\AS\\ 


«ec" ||S|| 


• Let v = v(x) be a function in a Hilbert space Ti, then 

IMI 2 = J v 2 (x)dx, 

• Let L be a differential operator L : H — > W where W is a Hilbert space, then 

ll^ll = su p nor- 

v(x)£H ll^ll 


Theorem [15] 

For any A 6 C mxn , the L 2 norm of A is given by 

\\A\\=<TM)^P h {A-A) 

where cr\ is the maximal singular value and p is the spectral radius. 


Corollary 

Using the above theorem and Parsevahs relation, the L 2 norm of a differential operator can be bounded 
by 

Ill'll < maxp^ ^L"(k)L(k)j, 

where L(k) denotes the symbol of the operator L. 

Appendix B. The Sensitivity Equation and the Adjoint Methods. Let L be a partial differential 
operator defining the following PDE: 


L{Q, &i, ■ * * , b\f) = 0, 


where Q denotes the state variables and bj denote the j’th design variable. Let I = I(Q,b\, • * • , 6m) be a 
functional. The sensitivity problem is to compute the derivatives 


(B-l) 



27 



A brief summary of the “sensitivity equation method” (also known as “direct- differentiation” ) and of the 
’’adjoint method” are given in the following. For a more extensive treatment of these methods see for example 
[16]-[17] for the sensitivity equation method and [12], [18] for the adjoint method. 

B.l. The Sensitivity Equation Method. The sensitivity equations are a linear set of M equations 
given by 


(B.2) 


L Q ^j + L b . = 0 for 1 < j < M, 


where the solutions ipj are the sensitivities of the state variables with respect to the design variables on all 
the domain: 


(B.3) 



In terms of the sensitivities, the derivatives of the cost functional with respect to the design variables are 
given by 


( B ‘ 4 ) 9j = oq'l’j + for 1 < j < M. 

B.2. The Adjoint Method. With the adjoint method the computation of the derivatives gj is done 
with the solution of only one linear PDE, rather than M, irrespective of the number of design variables. 
The adjoint equation, for the adjoint variable A, is a linear PDE is given by 

( B 5 ) L « A+ ^ =0 - 

where the operator Lq is the adjoint of Lq. In terms of the adjoint variable, A, the derivatives of the cost 
functional with respect to the design variables are given by 


(B.6) 


9j — ^6, A + 


dl 

dbj ' 


Appendix C. Lemma 1. 

Lemma 1: 

For a two states system (2.2), starting with errors Q\ and Q 2 in the state variables and ass umin g that 
the errors satisfy the quasi-linear approximation (2.3), a GGS sequential iteration (Ai,A 2 ) results in the 
following error relations: 


Qi +1 = fa** 

Q2 + 1 = fcQ^Q^QRl.Q^. 


Proof; 

The solution sequence of system (2.2) implies first solving R\ for Q\ and then solving R 2 for 

Q 2 as defined in Eqs.(2.4-2.5). The solution of Eq.(2.5) is given by 


(C.2) 


( Q 1 \ = -( R hQi 0 \ ( Rl \ 

\ Q 2 J \ ~ R 2,Q 2 ^ i Q 1 Ri Q 1 R 2 t Q 2 j \ R 2 J 
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Let (<2*> <25) denote the solution of the coupled system (2.2). The errors after a GGS iteration (2.4) are 
given by 


(C.3) 


Q" +1 = Ql-(Qi + Qi) = Qi-Qi 


QX 1 = Qi~ (<?2 + Q 2 ) = Qi- q 2 . 

Substitution of (C.3) in (2.3) gives 


(C.4) 


pn+l nn+1 

n kQi n hQ2 

pn+l pn+l 

^Qi /X 2,Q 2 


Q i - Q i 
<?2 — Q 2 


R n x +1 

i^ +1 


Substitution of (C.2) in (C.4) gives 

( R n X l ^ 


(C.5) 

where 

(C.6) 


K. = 


pn+l 

/t 2 


1 0 
0 1 


R 

K 


Tq\ R Vq\ ) ic(Q'\ 
Vq\ Kq\ ) \Q*J 


R 

R2 


1, Qi 0 \ / R\,Qi Rl, 

2, Qi #2,Q 2 / \ ^2,Qi #2, 


Q2 

Q2 


Using relations (2.3) and multiplication of Eq.(C.5) from the left by the matrix 

( pn+l pn+l \ * 

H 1,Qi ,Q a | 

pn+l pn+l / 

^Qi n 2,Q 2 / 

results in the following relation between the errors before and after a GGS iteration for MDA: 

( C 7) ( or 1 \ = IV 1 0 ) _ ( R ik 0 ^ ( R ^ Ri,Q, X] ( Qi 

V or 1 j Lv ° 1 / v -^k R ^ R -,k ^k ) v ^, Ql r 2 , q , )\ \ Q 2 


^Prom Eq.(C.7) it is concluded that 


(C.8) 


Qi +1 = -(Ri'Qj-'Ri'QtQi 


Q 2 ~ ^2,Q 2 ^'2 J Q 1 R 1,Qi RliQiQz* 

By the above two relations the following holds: 

(c.9) or 1 = -(i+,Q+i,Q 2 ) (wkx lR kk wkx'Rvk^Qr 1 

and 

(c.io) Qi = or 1 . 

Substituting of relation (C.IO) in relation (C.9) results in an expression for the error reduction in Q\: 

(c.ii) qx 1 = 
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Appendix D. Fourier Analysis of the Aeroelastic Problem. 

D.l. The Analysis Problem. The error equations are analyzed with Fourier analysis resulting in the 
symbols of the Jacobian operators, R%,Q r Then Parseval’s relation is applied to estimate the norm of the 
operators R%^Q j in the real space. 

The Flow Analysis Problem: R\(<j>,f3) = 0 


The analysis problem, R \ , is given by the following PDE (in this specific case the error equation is 
identical to the state equation): 


0 

II 

'O- 

<1 

in Q = {0 < x < 1 ; 0 < y < 1} 

II 

~0- 

on z = 0 

II 

0 

on x = 0 

0 

II 

<0 

on x — 1 

II 

0 

on z — 1 


where (3 = /?(x) is some non-constant L 2 [ 0, 1] function. 
The solution of Eq.(D.l) is given by 

(D.2) 4>(x,z)= 

k = — oo,k^0 

where 


(D.3) = — tanh(27r/c) cosh(27rA;z) + sinh(27rfcz). 

Since the function |t/>k(z)| has a maximum at z = 0 in the interval z e [0, 1], for all /c, we conclude that 
(D.4) |<^>(:r,z)| < |<^>(x,0)| 0<z<l 


and therefore, for the sake of error convergence, we will examine the error in <j) on the boundary only ( z = 0). 
By Eq.(D.2) on the boundary <f> = <j>(x) satisfies 


(D-5) 


#*, o) = -zr=- 

<t>z(x, o) = 52 


A tanh(27rfc) „j(2nkx) 
-oo,fc^0 Pk 2nk e 
°° A -i(2irfc*) 

/fc=-oo,fc#0^ e 


We conclude the following symbols of the Jacobians: 


(D.6) 


= -ssife 

Ri , a (&) = —2nik 


R h w(k) = —2nik 

RiMk) — 0. 


The Structure Analysis Problem: R, 2 (W,D) = 0 


The analysis problem, R 2 , is given by the following PDE (simply supported boundary conditions): 


(D.7) 


d xx {Dd xx W) = -p 
W{ 0) = W( 1) = 0 
d xx W{ 0) = d xx W{ 1) - 0. 
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Linearization of (D.7) results in the following error equation 


dxx{Dd xx W) + d xx {Dd xx W) = -p on 2 = 0 ; 0 < i < 1 

(D.8) 1^(0) = W( 1) = 0 

d xx W(0) = d xx W(l) = 0 . 

The errors W , D, and p are expanded by a Fourier series: 

00 00 00 

(D.9)W r (x) — W k e^ 2nk ^ ; D{x) = £ D k e^ 2 * kx) ; p(x) = £) p k e i{2 * kx) . 

k—~oc,k^0 k= — oo y k^Q k= —oc,k^Q 

Substitution of (D.9) into (D.8) results in an expression for the symbol of Rz- 

(D.10) R 2 = G w (D)W + G D (W)D+pt$ 

where 

G W (D) = D(2nk) 4 - 2i{d x D)(2Trkf - ( d xx D){2irk ) 2 
(Dll) G d (W) = -(2nk) 2 d xx W + 2(2mk)d xxx W + d xxxx W 

p $ — 2nik. 

Therefore, the symbols of the Jacobians are give by 

R 2 , 4 >{k) = 2irik ; R 2> w(k) = Gw(D) 

i*2,a(*) = 0 ; R 2 M^) = Gd(W). 

D.2. The Adjoint Problem. 

The Adjoint Flow Problem 

The symbol of the adjoint flow operator P\ on the boundary is given by (Pi (A, 77, <f) , W) — d z X - d x rj — 
27i^xx0-73 d x W) 

(D.13) Am = - - <*r*>(« + 73*) + 27,(2,*,)^. 

therefore the symbols of the Jacobians are given by 

A, </>(£) = 2 ji( 2 nk ) 2 Pi,w(k) = -73(2 nik) 

( D14 ) A.a (*) = - 5s|fe A,„(*) = -2*ifc 

A,«W = 0 A,d (*0 = 0. 

The Adjoint Structure Problem 

The symbol of the adjoint structure operator P 2 on the boundary is given by (P, (A, f], <p, D) = d xx (Dd xx T])+ 

d x \ + ^d x 4>) 

(D.15) P 2 (k) = G w {D)f, + G d { V )D + (27rik)(X + 73^), 

therefore the symbols of the Jacobians are given by 

A,tf.(A = 73(27 rik) P 2 ,w(k) = 0 

(D.16) p2,\(k) = 2 t rik A ,„(*) = G W (D) 

P2,a(k) = 0 A.d(*) = G d (»7). 
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[Ai,A<2 



n 


Fig. 6. Actual convergence rate versus theoretical upper bound for the GGS (sequential) solution of the non-linear algebraic 
system. 



Fig. 7. Actual convergence rate versus theoretical upper bound for the MD feasible (left) and single discipline feasible 
(right) sensitivity Solution of the non-linear algebraic system . 
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Fig. 8. Actual convergence rate versus theoretical upper bound for the MD feasible (left) and single discipline feasible 
(right) optimization solution of the non-linear algebraic system. 
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Fig. 9. Actual convergence rate versus theoretical upper bound for the GGS (sequential) solution of the aeroelastic problem. 
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Fig. 10. Actual convergence rate versus theoretical upper bound for the MD feasible (left) and single discipline feasible 
(right) sensitivity solution of the aeroelastic problem . 
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Fig. 11 . Actual convergence rate versus theoretical upper bound for the single discipline feasible optimization solution of 
the aeroelastic problem. 
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